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A robust all-at-once multigrid method for the Stokes 
control problem 


Stefan Takacs 


Abstract In this paper we present an all-at-once multigrid method for a 
distributed Stokes control problem (velocity tracking problem). For solving 
such a problem, we use the fact that the solution is characterized by the 
optimality system (Karush-Kuhn-Tucker-system). The discretized optimality 
system is a large-scale linear system whose condition number depends on the 
grid size and on the choice of the regularization parameter forming a part of the 
problem. Recently, block-diagonal preconditioners have been proposed, which 
allow to solve the problem using a Krylov space method with convergence rates 
that are robust in both, the grid size and the regularization parameter or cost 
parameter. In the present paper, we develop an all-at-once multigrid method 
for a Stokes control problem and show robust convergence, more precisely, 
we show that the method converges with rates which are bounded away from 
one by a constant which is independent of the grid size and the choice of the 
regularization or cost parameter. 

Keywords PDE-constrained optimization • all-at-once multigrid methods • 
Stokes control 


1 Introduction 

In the present paper, we consider the following model problem (distributed 
velocity tracking problem). Let 17 C be a bounded domain with d G {2, 3}. 
Find a velocity field u G [i7^(l7)]‘^, a pressure distribution p € and a 
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control (force field) / G such that the tracking functional 

1 rv 

is minimized subject to the Stokes equations 

—Au + Vp = / in 17, 

V • M = 0 in 17, 

M = 0 on dQ. 

The cost parameter or regularization parameter a > 0 and the desired state 
(desired velocity field) ud G [L'^{fl)Y are assumed to be given. To enforce 
uniqueness of the solution, we additionally require f^p dx = 0. 

Here and in what follows, Lp‘{fi) and denote the standard Lebesgue 

and Sobolev spaces with associated standard norms |j • \\l'^(^q) and || • Hz/qr?), 
respectively. 

The main goal of this work is to construct and to analyze numerical meth¬ 
ods that produce an approximate solution to the optimization problem, where 
the computational complexity can be bounded by the number of unknowns 
times a constant which is independent of the grid level and the choice of the 
parameter a, in particular for small values of a. 

The solution of the optimization problem is characterized by the Karush- 
Kuhn-Tucker-system (KKT-system). As we are interested in good approxima¬ 
tions of the solution, the discretization of the KKT-system leads to a large-scale 
linear system. This linear system will be solved with multigrid methods be¬ 
cause they are one of the fastest known methods for such problems. Originally, 
multigrid methods have been designed and analyzed for elliptic problems. They 
also work well for saddle point problems (like KKT-systems) and have gained 
growing interest in this area, see, e.g., H] and the references cited there. 

The unknowns of the discretized KKT-system for a PDE-constrained opti¬ 
mization problem can be partitioned into primal variables and dual variables. 
In our case, the primal variables consist of state variables (the velocity field u 
and the pressure distribution p) and control variables (the force field /). The 
dual variables are the Lagrange multipliers which are introduced to incorpo¬ 
rate the constraints. One approach to solve such problems is to apply multigrid 
methods in every step of an overall block-structured iterative method to equa¬ 
tions in just one of these blocks of variables. Such methods have been proposed, 
e.g., in [UlinillT]. 

Another approach, which we will follow here, is to apply the multigrid 
idea directly to the (reduced or not reduced) KKT-system, which is called an 
all-at-once approach. Such methods have been proposed and discussed for the 
elliptic optimal control problem, e.g., in nang. 

In this paper we present a convergence proof for multigrid methods based 
on the classical splitting of the analysis into smoothing property and approx¬ 
imation property, see [7j. 
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This paper is organized as follows. In Section]^ we will introduce the op¬ 
timality system and discuss its discretization. In Section we will introduce 
an all-at-once multigrid approach. Its convergence will be proven in Section]^ 
Numerical results which illustrate the convergence result will be presented in 
Section [sl In Section |6] we will close with conclusions. 


2 Optimality system and discretization 

For setting up the optimality system we need the space i?g(I?), the space of 
functions in vanishing on the boundary. Moreover, we need the space 

Lg(l7), which is the space of functions in with mean value 0, i.e., 

Lg(l7) := In e : J v = 0 

Both spaces are equipped with standard norms, i.e., || • •= II ’ llf/qfi) 

and II • ||l 2(^3) := || • ||i2(^2). 

The solution of the problem is characterized by the Karush Kuhn Tucker 
system (KKT-system), which reads as follows, cf. [18] and others. 

Find e [H^{f2)]‘‘xLl{n)x[L^{f2)]'^x[H^{f2)]‘^xLl{f2) such 

that 


{u, u)L^(n) -I- (VA, X/u)l 2 i^q) -I- (^, V • u)L^(n) = {ud, 

(V • A,p)i2(j7) = 0 

«(/) f)L^(n) - /) l 2 ( j 7 ) = 0 

(Vm, VA)l2(i7) -I- (p, V • A)i2(j7) - (/, A)i2(j7) = 0 

(V • U, p)l2(j7) = 0 

holds for all {u,pJ,X,ll) G x Ll{n) x [L^{f2)]‘^ x [H^{n)]‘^ x L^n). 

The third line of the KKT-system directly implies / = a“^A. 

This allows to eliminate the control /, which leads to the reduced KKT- 
system, which reads as follows. Find x := {u,p,X,p.) £ X := [Fig (17)]'^ x 
Lg(J7) X X Tg(l7) such that 

+ (VA, Vu)i2(j7) -b (p, V • u)L^n) = {ud,u)l^o) 

(V • A,p)i2(i7) = 0 
(Vu, VA)L2(f2) -t- (p, V • X)L2(n) — a“^(A, A)l2(i7) = 0 

(V • u,jj.)L2(a) = 0 


holds for all (u,p, A,/i) G X. 

Certainly, the KKT-system can be rewritten as one variational equation as 
follows. Find x G X such that 


B{x,x) = J-{x) for all x G X, 


(1) 



4 


Stefan Takacs 


where 


{u,p,X,fl)) := {u,u)l2(s7) + (VA, Vu)l2(^3) + (/x, V • u)l^(o) 

+ (V • A,p)l2(i 7) + (Vu, VA)L2(f2) + (p, V • A)L2(f2) — a~^(A, A)l2(i7) 

+ (V • u, A)i 2 (x 2 ) and 
T{u,p,X,il) := (uD,u)L2(ay 

We are interested in finding an approximative solution for equation Q. 
Both, the proposed solution strategy and the convergence analysis, follow the 
abstract framework introduced in [16j . The conditions, (Al), (Ala), (A3) 
and (A4), mentioned in the present paper are the same conditions as in [16j . 
For simplicity, we introduce the following notation. 

Notation 1 Throughout this paper, C > 0 is a generic constant, independent 
of the grid level k and the choice of the parameter a. For any scalars a and b, 
we write a ^ b (or b ^ a) if there is a constant C > 0 such that a < C b. We 
write a ^ b if a b < a. 


The following property guarantees existence and uniqueness of the solution. 


(Al) The relation 


holds for all a; G A. 


Ikllx < sup 

O^xGX 


B{x, x) 


<lkllx 


In [TB] it was shown that condition (Al) is satisfied for X :=Y xY, where 
Y := U X P, U := [i?Q(l7)]‘^, P := Lq( 12), equipped with norms 

llxlll := ||(«,p,A,/i)||| := ||(u,p)|||. +a-i||(A,M)|||-, 


where 


iu,p)ry:= \\u\\fj + \\p\\l, 


W\u ■— \W\L^(n) 
Ibllp := sup 


and 


a{p,X 




d w 


L2(r?) 


a 


1/21 


W\ 




Using the following notation, we can express the norms in a nicer way. 

Notation 2 For any Hilbert space A, the symbol A* denotes its dual space 
equipped with the dual norm 


|m||a* := sup 


(m, w) 


o^weA II^IIa 


where {u, •) := u(w) denotes the duality pairing. 

For any Hilbert space A and any scalar o > 0, the symbol a A denotes the 
space on the underlying set of the Hilbert space A equipped with the norm 
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For any two Hilbert spaces A and B, the symbol An B denotes the space on 
the intersection of the underlying sets, {u G An B}, equipped with the norm 


Ihll 


2 

AnB 


IhllUlkll 


2 


and the symbol A + B denotes the space on the algebraic sum of the underlying 
sets, {mi + U 2 ■ ui G A,U 2 G B}, equipped with the norm 


\A+B ■- 


inf 

UiGA,U2GB,U — Ui-\-U2 


\uin+\\u2\\i. 


The spaces A*, a A, AnB and A + B are Hilbert spaces. The fact that A* 
is a Hilbert space follows directly from the Riesz representation theorem, see, 
e.g.. Theorem 1.2 in [1]. The fact that oH is a Hilbert space is obvious and 
for the latter two see, e.g.. Lemma 2.3.1 in [5]. 

We immediately see, that the norm on U can be rewritten as follows 


Iklltf — \\u\\L2(n)na^/2H^{n)- 

To reformulate the norm || • ||p, we need a regularity assumption. 

(R) Regularity of the generalized Stokes problem. Let / G and g G 

HQ(f2) n L§(12) be arbitrarily but fixed and (u,p) G [iLo(l7)]‘^ x LQ(f2) be 
the solution of the Stokes problem, i.e., such that 

(Vu,Vu)/,2(f2) + (p, V • M)p2(f2) = (/, M)L2(i2) 

(V ■ U,p)^2(f2) ={9,P)l2(0) 

holds for all {u,p) € x Lq( 12 ). 

Then {u,p) G x H^{Q) and 

+ Ib||p-I(f2) < Il/lli2(j7) + Il5llffi(r2)- 

This condition is satisfied for convex polygonal domains, see Lemma 2.1 
in [13] which is a direct consequence of Theorem 2 in |S] . 


Lemma 1 If (R) is satisfied, then 

Wp^MI aH^{n)+a^/^L^(S7) 

holds for all p £ Lq(12). 

This lemma was shown in Theorem 3.2 in HD under a regularity assump¬ 
tion, which is weaker than regularity assumption (R). 

The discretization of problem 0 is done using standard finite element 
techniques. We assume to have a sequence of girds obtained by uniform re¬ 
finement. On each grid level k, we discretize the problem using the Galerkin 
approach, i.e., we have finite dimensional spaces Xk C X and consider the 
following problem. Find Xk G Xk such that 


B{xk,Xk) = X{xk) for all Xfc e Xfc. 


(2) 



6 


Stefan Takacs 


Using a nodal basis, we can represent this problem in matrix-vector notation 
as follows: 

•^k = f (3) 

Here and in what follows, any underlined quantity, like is the representation 
of the corresponding non-underlined quantity, here with respect to a nodal 
basis of the corresponding Hilbert space, here Xk- 

Existence and uniqueness of the discretized problem is guaranteed by the 
following condition. 

(Ala) The relation 


\\xk\\x < sup 

O^XkGXk 


l3{xk,Xk) 


^ W^kWx 


holds for all Xk G X^. 

Due to the fact that the model problem is indefinite, condition (Al) does not 
imply condition (Ala). For the Stokes problem itself, it is well-known that 
such a condition (also known as discrete inf-sup eondition) can only be guar¬ 
anteed if the discretization is chosen appropriately. The same is true for the 
Stokes control problem. Fortunately, we can show the discrete inf-sup condi¬ 
tion (Ala) for the Stokes control problem based on pre-existing knowledge 
on the discrete inf-sup condition for the Stokes problem. This allows to show 
that all discretizations which are suitable for the Stokes flow problem are also 
suitable for the Stokes control problem. 

We choose the space X^ as follows: 


Xk ■■= Yk X Yfe where := Uk x P^. 


and the choice otUk Q U = [i7g (12)]'^ and C P = Lg)!?) is discussed below. 
Note that Xk has product structure and that the state and the adjoined state 
(Lagrange multipliers) are discretized the same way. The same has already 
been done for optimal control problems with elliptic state equation, cf. [16] 
and many others. 

Due to the fact that the grids are obtained by uniform refinement, the 
discrete subsets are nested, i.e., Uk Q Uk+i and Pk C Pk+i- Therefore, also 
Xk C Xk+i holds. 

The next step is to show condition (Ala). We have seen that the analysis 
done in [TB], applied to the infinite dimensional spaces, shows condition (Al). 
If the analysis done in [TB] is applied to the discretized spaces, we obtain that 


\\xk\\Xk 


< sup 


B{xk,Xk) 

WXkWxk 


^ W^kWxk 


(4) 


is satisfied for all Xk & Xk, where 

ll^fellxfc ■ \\{'^k,Pk, ^k, ikk)\\xk ■ II Pfe) II Y), Y CX IK-^fc) Mfe) || Yj, j 


\\{uk,Pk)rn ■■= Ml + WPkWl,, 


Ibfellp, := sup 


a(pfc,V • Wfc)|2(j7) 


O^WkGUk Il'*^/c|b2(^) +' 




and 
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II • 11^ is as above. 

Note that this is not condition (Ala), as the norms || • ||p and || • ||pj, are 
not equal. For showing condition (Ala), it suffices to show that these two 
norms are equivalent which implies also the equivalence of the norms || • ||x 
and II • \\xk - This can be shown using the following condition. 

(S) The discretization of P is iJ^-conforming, i.e., Pk Q and the weak 

inf-sup condition 


(V • Uk,Pk)L^(n) 

-n—n- 

O^UkGUk \\Uk\\L^(n) 


^ l|Vpfe||L2(j7) 


holds for all Pk & Pk- 


Lemma 2 Assume that the discretization satisfies condition (S). Then con¬ 
dition (Ala) is satisfied for the model problem. 

Proof Lemma 2.2 in m states (provided that (S) is satisfied) that || • ||p ~ 
II • ||p^ is satisfied. A direct consequence is || • ||jc ~ || • \\xu- Therefore, condi¬ 
tion Q implies condition (Ala). □ 

Note that condition (S) is a standard condition which ensures that the 
chosen discretization is stable for the Stokes problem. In HE] it was shown 
that condition (S) is satisfied for the Taylor-Hood element (PI — P2-element) 
for polygonal domains where at least one vertex of each element is located in 
the interior of the domain. Here and in what follows we assume that the prob¬ 
lem is discretized with the Taylor-Hood element and that the mesh satisfies 
the named condition. 


3 An all-at-once multigrid method 

The problem shall be solved with an all-at-once multigrid method. The ab¬ 
stract algorithm for solving the discretized equation @ on grid level k reads as 
follows. Starting from an initial approximation , one iterate of the multigrid 
method is given by the following two steps: 

— Smoothing procedure: Compute 

:= - Ak for m = 1,..., 12 

with ■ The choice of the smoother (or, in other words, of the 

preconditioning matrix A'jl^) will be discussed below. 

— Coarse-grid correction: 

- Compute the defect := f ^ — Ak i>±k’'''^ restrict it to grid level 
k — 1 using an restriction matrix Ik~^- 
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— Solve the coarse-grid problem 




^k-l 


( 5 ) 


approximatively. 

— Prolongate ^ to the grid level k using an prolongation matrix /^_i 
and add the result to the previous iterate: 


™(1) 


+ Ik-1 P 


( 1 ) 

k-1' 


As we have assumed to have nested spaces, the intergrid-transfer matrices 
and can be chosen in a canonical way: I^_i is the canonical embedding 
and the restriction is its transpose. 

If the problem on the coarser grid is solved exactly (two-grid method), the 
coarse-grid correction is given by 

¥-k^ ■= 4°’"^ + ll~^ (4 - A ■ ( 6 ) 

In practice the problem © is approximatively solved by applying one step 
(V-cycle) or two steps (W-cycle) of the multigrid method, recursively. On the 
coarsest grid level (fc = 0) the problem ([^ is solved exactly. 

To construct a multigrid convergence result based on Hackbusch’s splitting 
of the analysis into smoothing property and approximation property, we have 
to introduce an appropriate framework. 

Convergence is shown on the spaces , which are equipped with an L^-like 
norms ||| • |||o,fc, which are defined a follows: 

llla^fclllo.fc := kfcllL ■= {IIkXk,Xk)p, 

where 


/ Pa,kMu,k \ 

Pa,kI^P,k 

a~^(pa^kMu,k 


( 7 ) 


and ipa,k := 1 -b and Mjj k and Mpk are the mass matrices, repre¬ 

senting the L^-inner product on Uk and Pk, respectively. Based on the norm 
III • lllo.fc, we can introduce the residual norm ||| • ||| 2 ,/c using 


llkfcllb.fc 


sup 

Xk^Xh 


B{xk,xk) 

|||ife |||o,fe 


Smoothing property and approximation property read as follows. 

— Smoothing property: 

Ill4°’''^ - 4ll|2.fe < Tl{p)\\\Xk^ - 4ll|o.fc (8) 

should hold for some function ri{v) with lim,y_>oo vip) = 0 - Here and in 
what follows, a:^ € Xk is the exact solution of the discretized problem ©. 
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— Approximation property: 


4^^ -4lllo,fc < Ca||| 4°’‘'^ - a; 


fcll|2,fe 


(9) 


should hold for some constant Ca > 0. 

It is easy to see that, if we combine both conditions, we obtain 

ll|a;4 - 4lllo,fc < - 41110,fc, 


where q{iy) = CArjiv), i.e., that the two-grid method converges for v large 
enough. The convergence of the W-cycle multigrid method can be shown under 
mild assumptions, see e.g. [7]. 

The choice of an appropriate smoother is a key issue in constructing such 
a multigrid method. Here, we introduce one smoother which is appropriate 
for a large class of problems including the model problem: normal equation 
smoothers^ cf. [S], which read as follows. 


(0,m) 




:= X 


(0,m-l) 




T Ck ^AkC-f. 


lilk 


-A 


k X 


(0,m — 


^k 


■’) 




for m = 1,..., i/. 


Here, a fixed r > 0 has to be chosen such that the spectral radius p{TA^^Ak) is 
bounded away from 2 on all grid levels k and for all choices of the parameter a. 
Using a standard inverse inequality, one can show that 


\Xk\\x 


< 


\^k llo.fc 


is satisfied for all Xk G Xk- Based on this result, using an eigenvalue analysis 
one can show the following lemma, cf. [S]. 

Lemma 3 The damping parameter t > 0 can be chosen independent of grid 
level k and the choice of the parameter a such that 

T p{Af^Ak) < 2 - e < 2, 

holds for some constant e > 0. For this choice of t, there is a constant Cs > 0, 
independent of the grid level k and the choice of the parameter a, such that 
the smoothing property Q is satisfied with rate 

V{^) ■= 

Certainly, the iteration procedure Q should be eSicient-to-apply. Using 
the fact, that the mass matrices and Mp^k in Q and their diagonals are 
spectrally equivalent under weak assumptions, for the practical realization of 
the smoother these matrices can be replaced by their diagonals. 
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4 A convergence proof 


The proof of the approximation property is done using the approximation the¬ 
orem introduced in m which requires besides the conditions (Al) and (Ala) 
two more conditions (conditions (A3) and (A4)) involving, besides the Hilbert 
space X, two more Hilbert spaces ■= (-^-,|| • l|jf-,fc) and X^^k '■= 

(A_|_, |j • IIjc,,. fc), which are chosen as follows. 

As weaker space, we choose A_ := IT x F_, where Y- := U- x P_, 
U- := [L'^{[2)]‘^ and P_ := [i7g(l7) n Lq( 17)]*, equipped with norms 

lla^llL., := II(«,P,A,m)||L,. := ll(«,P)llL . + «-'I1 (A,m)IIv_.., 

ll(«,P)llv_,. := lhllL.. + lbllp_.., 

Ibll(7_,fc lbll[pi(r2)+a-i/2_L2(^2)]* and 

lbllp-,fc ■= ^fc^lbllfa-iLg(r2)na-i/2pi(l7)]* 

Note that dual spaces are (A_)* := (IT)* x (H_)*, where (IT)* = (H-)* x 
(P_)*, {U-Y = [L‘^{n)Y and (P-)* = H Lq( 17), equipped with norms 

ll-^ll(x_,.). := ll(/,<?,C,x)ll^x_..)* ll(/-5)ll(r_..). +«ll(C,x)ll^r_..)., 

\\f\\fu-,k)* ^fell/llpi(r 2 )+a-i/ 2 i 2 (^ 2 ) and 

ll5ll(P_^^)* := ^k\\9\\a-iLl(n)na-^/2H^(n)- 

As stronger space, we choose A+ := 1+ x F+, where F+ := P+ x P+, 
U+ := [P^(17) n HQ{f2)Y and P+ := n Po(l7), equipped with norms 

lblll++ := ll(«,P,A,M)||i,,, := \\{u,p)fr^,,+a-Yi\,9)\\h.^ 
ll(«,P)lly++ := lbllp++ + lbllp++, 

lbllp++ ■= ^fclbllpi(r 2 )nai/ 2 // 2 (f 2 ) and 
lbllp++ ■= ^fclbllaP2(f2)+Q,i/2//l(r3)- 
The additional conditions read as follows. 


(A3) On all grid levels k, the approximation error result 


inf lb - ccfellx < |b||x+fc for all a; € A+ 

xk&Xk 


is satisfied. 

(A4) For all grid levels fc, all P € (A-)* the solution xjr g A of the problem. 


find X G X such that B{x,x) = P(i) for all x G X, (10) 


satisfies xjr g A+ and the inequality 


lbp||x+,fc < l!A|l(v:_g^).. 


(11) 
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Based on these assumptions, the following theorem shows the approxima¬ 
tion property. 

Theorem 1 Let for k = 0,1,2,... the symmetric matrices Ak be obtained 
by discretizing problem © using a sequence of finite-dimensional nested sub¬ 
spaces Xk-i C Xk C X. Assume that there are Hilbert spaces X+ C X C X_ 
with mesh-dependent norms || • |jjf+ t? II ■ II a and || • ||js:_ ^ such that the condi¬ 
tions (Al), (Ala), (A3) and (A4) are satisfied. Then the coarse-grid cor¬ 
rection § satisfies the approximation property 

Il4 -4 IIa_., < sup ^ II II - (12) 


where the constant Ca only depends on the constants that appear (implicitly) 
in the named conditions. 

For a proof, see [16], Theorem 4.1. 


Theorem 2 Condition (A3) is satisfied. 

Proof This proof is analogous to the proof of Theorem 4.2 in El- However, to 
keep this paper as self-contained as possible, we give a proof of this theorem. 

Note that it suffices to show approximation error results for the individual 
variables separately. Using a standard interpolation operator Hk : [H^(l7)]'^ — 
Uk, we obtain for the velocity field u 

||w — 7Tfcu||^2(f2) ^fe||u||ffi(j7) \\n — T[ku\\'jji(^Q^ < hl\\u\\‘)j 2 i^Qp 


for all u € and therefore 

inf ||w - Ufcll^ < ||u - 77fcM||^ = ||u - HkuWl^f^a) + a^/^||w - 

'^k k 

^ ^k {\\n\\'Hi(o) + a^'^^l|w|lif2(^)^ = IlMlI^.k- 


The same can be done for the adjoined velocity A. Also for the pressure dis¬ 
tribution p we can do a similar estimate. The estimates 


inf \\p-pk\\h(^n) 
PkePk ' ' 


< 


^felbll 111(11) 


and 


inf 

Pk&Pk 


\P-Pk\\m{o) 


^ <hi\\pr 


H^(n) 




12 


Stefan Takacs 


are standard approximation error results which imply 

lb - Pfellp = lb " Pk\\lm{n)+c.U2mn) 

= \\<ll\\lH^a) + \\<l2\\ai/2L2(0) 

qiGH^n) 
q2^L‘^ {O) 
qi+q 2 =p-Pk 

= IIpi + inf lb2-P2,fc||^l/2i2(^2) 

pi^H^(n)pi,kePk ' ' p2,kePk ' 

P2eL^{n) 

Pl+P2=P 

- inf Ibi-Pi,fcllapi(r3) + inf lb2-P2,fc||pi/2i2(^j) 

pi^H^(n)pi,kePk ' ' p2,kePk '■ ’ 

P2GH\n) 

Pl+P2=P 

^ eH^(O) 11^211^1/2^1(17) = ^felbllap 2 (i 7 )_,_„i/ 277 i(i 7 )- 

P2effHf2) 

Pl+P2=P 

The same can be done for the adjoined pressure /i. This finishes the proof. □ 

For showing (A4), we recall Theorem 4.6 in [T3j on the regularity of the 
generalized Stokes problem. For this purpose, we need a regularity assumption 
for the Poisson problem with homogeneous Neumann boundary conditions. 

(Rl) Regularity of the Poisson problem. Let g G and p G H 

Lq{Q) be such that 

(Vp, Vp)H^n) = ig,p)L' 2 {n) for all p G H^{[2) n Lo(f^)- 
Then p G and Ibl|p 2 (r 2 ) < |b||L 2 (i 7 ). 

Such a regularity assumption can be guaranteed for convex polygonal domains 
(see, e.g., 0). 

Theorem 4.6 in m directly implies the following theorem. 

Theorem 3 Suppose that the regularity assumptions (R) and (Rl) are sat¬ 
isfied. Let f G [L‘^{n)]‘^ and g G iLg (17) H Lq(17). The solution of the problem, 
find {u,p) G Y such that 

a"^/^(M, m)l 2(^3) + (Vu, V{t)i2(^3) + (p, V • u)z,2(f2) = (/, u)i2(^3) 

(V • U,p)i 2 (f 2 ) = {g,p)L 2 (n) 

for all {u,p) G V, satisfies {u,p) G 1+ and the inequality 

lblla-i/27/i(r2)nP2(r2) + lbllQi/2//2(^)_|_77i(^) 

^ ll/llai/2i/i(f7)+_L2(f2) + Iblla-1/2L2(12)P,771(^) 

is satisfied. 

Proof We choose the parameter /3 (which occurs in [T3]) to be fi := . □ 
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Lemma 4 Suppose that assumptions (R) and (Rl) are satisfied. Let T S 
(X_)* he arbitrarily but fixed. Then, xjr, the solution of (10), satisfies xjr S 
X_|_ and the bound 


Proof Let P{u,p,X,il) := (/,w)L 2 (f 2 ) + {g,p)L'^(n) + (C,X)L^(f7) + (x,fi)L^o), 
where /,C G [^^(12)]'^ and g,x& n Lg(12). 

Let f ■= f — ujr + and C := C + ct“^Ajr + a~^l'^ujr. Then we can 

rewrite the KKT-system as follows: 

(VAjf, Vu)i2(j^) + a + (/tjr, V • 'i()z,2(^) = (/, u)L2(r2) 

(V • Ajr,i5)z,2(f2) = {g,p)L^(n) 


and 


(Vmjt, VA)i 2 (f 2 ) + a A)i 2 (i 7 ) + (pjt, V • X)L 2 (n) = (C A)L 2 (f 2 ) 

(V • Mjr,/1)L2(17) = {x,P)L2(a)- 

As / G [L2(12)]<i, g e Ri(12)nL2(r2), ^ g [L2(12)]'^ and x G Hl{n)nLl{L2), we 
obtain using Theorem]^ that xjr G X+ and the following bounds are satisfied: 

^ W f — + Oi ! ^J^\\oL^/^Hl(n)+L'^{Q) \\9\\a-^/^Ll(Q)r]Hl{Q) 

and 

l|■a^llQ-l/2_f/l(^2)n_f^2(^2) + lb-^llai/2ii-2(^)_|_^l(^) 

<||C —a ^Ajr + a ^^^uj='||q1/2//1(^2)_|_j;,2(i7) + ||xIIq-i/2l2(^)p,^i(^). 

We can combine these two estimates and obtain 

+ lb-^llc<_fl2(J7)+aV2/i-i(j7) 

^ ll/llffi(12)+c<-i/2L2(^2) + ll5lL-iLg(12)a-i/2nili(f3) 

+ '^IICII//i(i7)+Q,-i/2i2(f2) + '^llxlU-iL2(r2)na-i/2//i(f2) 

+ lb.7='llffi(t2)+a-i/2L2(f2) + <a \\^y^\\Hfin)+oi-^/^L‘^{n)- 

Note that ||ujF||^i(^ 3 )+„-i/ 2 i 2 (^ 3 ) < ||ujr||^i(^ 2 ) = \\ujr\\mia) holds because of 
ujr g As the analogous holds also for Xjr, this finishes the proof. □ 
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To show condition (A4), we have to bound 
from above. For showing such a result, we need some notation. 

As i?g(l7) is dense in L^(l7), for u S the function —Au G 

[L^(17)]'^ can be approximated by some function w'^ G [HQ{n)]^ such that 

li - < e. 

So, we can introduce an operator —A'^ : —>■ [i?g(l7)]‘^ such that 

II — Au — (—A'^)m||| 2 (^ 2 ) ^ £• 


Analogously, we introduce the operator V' : H^{n) —)• [iJQ(l7)]‘^ such that 

Lemma 5 Let T G (A_)* and let xjr = (ujr,pjr,\jr^ he the solution 
of (10). Then xjr satisfies the estimate 

< ll-^ll(x_,.).l|a:^IU,,.. (14) 


hi 


Proof Let T{u,p,X,fi) := if,u)L 2 (^n) + {g,p)L^(n) + {C,^)L2{n) + 
where /,C G [L^(n)]‘^ and g,x G H^{^) n Lo(n). 

The idea of this proof is to show that for all e > 0 there is some x’^ G X 
such that 


X{x'') - B{Xjr,X^) 

^ ^fe^ll-^ll(x_,fc)*||a;x||x+,fc - l|Mx||ffi(r2) - a“^ll V||ffi(r2) (15) 
+ e(ai/^ + a-^/'^)hl^{\\xjr\\x^^^ + ||^||(x_,fc)*) + e^- 

Note that the left-hand-side of the inequality is 0. Therefore, this would be suf¬ 
ficient to show the statement of the lemma, as e > 0 can be chosen arbitrarily 
small. 

In the following, we show that (EH) is satisfied for the choice x'^ := {—A'^ujr, 
-X-X^pjr, A'^Xjr, X Wb estimate the individual summands of X{x'^) — 

B{xx,x^) separately. For the first one, we obtain 

- (mJF, -A'^Mjr)i 2 (^ 2 ) < -{ujr, - Aujr) + £ 11 11 ^2 

= -(yujr, Vwjr)i2(j7) -I- e||Mjr||i2(i7) < “ 11 Wx 11 1 (r?) + elifc ^Ikxllx+.fc 

due to the fact that ujr g [i?^(l7) n i/g (17)]'^ and due to Friedrichs’ inequality. 
The same can be done for a“^(Ax, A'^Ax)l 2 (i 7 ). 

For the next two summands, 

— (Vux, VA'’Ax)L2(f2) — (VAx, V(—A')'ux)l2(i7) 

= (Aux, A'^Ax)l2(j 7) — (AAx, A'^wx)l2(j7) 

< (Aux, AAx)l 2(^) - (AAx, Aujr)L2(Q) + e(|l Aux||L2(f2) -t- 1 | AAx||L2(f2)) 

< e(l|ux||ff2(^2) -I- ||Ax||ji'2(f2)) < ehfc ^(a^A -f a“^^'‘)||a;x||x+,fc 
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is satisfied due to the fact that A'^ maps into 
For the next two summands, we obtain 

— (V • ujr, V • V’^/rjr)L2(r2) — (V • {—A^)uj^, 

= (VV • ujr, - {A^uj^, 

< (VV • Ujr, V’^^jr)^2(j7) — (A’^Ujr, V‘^/Xjr) ^2 ( j7) + £ 11 11 ^2 ( 

< —(VUJF, VV'^//jr)j;,2 — {Aujr, fijr) ^^2 + e{\\Aujr\\ ^^2 + ||V'^Pjr||^2 + c) 

= —(VUJF, VV'^/rjr)j;,2 + (VUJF, VV'^/Xjr)^2 + e{\\Aujr\\^^2 + II V'^pjr||^2 + c) 

< e(l|wjF||i/2(f2) + Wpj^Wh^o) + 2 e) < + a“^'^^)||a;jr||x+,fc + e^- 

The same can be done for -(V • Ajr, -V • V^px)L2(fi) - (V • A'^Xjr,pjr)i^ 2 (a). 
Let /2 G and fi := / - /a. Then 

(/l; -^^'ajr)L2(i^) < ||/l||a-i/2L2(r2)||Mjp||ai/2/i-2(j7) + ea^/'‘||/i ||„-i/2i2(f2) 

holds as well as 

(/2) —^'^Ujf)l2(j7) < (V/ 2 , Vujf)l2(i 7) + e||/2||L2(i7) 

^ ll/2||ffi(r2)||ujF||ffi(r2) + e||/2||Hi(i7)- 

This implies 

(/, 

+ ehf^\l + a^/‘^)\\T\\ 

(X-A*- 

Let p 2 G H‘^{Q) and pi := pj^ — P 2 G iL^(l7). We have 

{g, —V • V'^pi)i 2 (j^) = (Vg, y^pi)L^(o) 

^ ll5llo!-i/2/i-i(j7)|bi||Qi/2^i(j7) + e||5||_f/i(r2)- 

Moreover, using g G iLg(l7), we have also 

{9, -V • V'^p2)L2(r2) = (V5, Vp2)L^(n) < (Vg, ^P2)l^(q) + e||g|l/i-i(r2) 

= -(5i V • Vp2)L2(r2) + e||5llffi(^2) < ll5lla-iL2(f2)|b2|U_f/2(f2) + e\\g\\H^(a) 
and therefore 


{9, -V • V'^pjF)L2(r2) 

^ ll5lla-iL2(i7)na-i/2ii-i(r2)|bxl|ai/2/i-i(r2)+a_f/2(i7) + ellffll 
^ \\9\\a-^L^(n)na-^/^H^{0)\\Px\\a^/2H^{0)+aH'^(n) + ll-^ll (X_ _fc)* 


is satisfied. The same can be done for {(^, A"^Xjr) 1^2 and (x, V • 

Combining these results, we immediately obtain (15), which finishes the 
proof. □ 
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Theorem 4 Condition (A4) is satisfied. 


Proof By combining (|T3|) and (14), we obtain 


for some constant C > 0 (independent of k and f) which implies 

||xx||x+,. < i [C+V^C + C^) ||^||(x_..)*, 

i.e. which finishes the proof. □ 

So, we have shown condition (A4). So, Theorem implies the approxi¬ 
mation property. Note, that we have now shown the approximation property 
in the norm || • ||x_ i.e., (12). The next step is to show the approximation 

property in the norm-pair ||| • |||o,fc and ||| • ||| 2 ,fc, he., (§. 

To show (§, the following lemma is sufficient. 

Lemma 6 The inequality 


Ilkfclllo.A 


< 


\Xk\\x_ 


(16) 


is satisfied for all G A^. 


Proof The proof of this lemma is based on Lemma 4.7 in [T3]. Lemma 4.7 
states (in the notation of the present paper and for the choice j3 := 
that 

-1/2| m ,_|||2 ^ „-1/2| L ,_ ||2 


a 




is satisfied for all yk &Yk- As for Xk = {ykjtpk) & Xk = Yk x Yk both. 


lll(j/fc.V'fc)lllo.fc = lllyfclllr.o.fc+« ilV'/clllr.o.fc 


and 


ll(yfe,'0fc)ll 


2 






is satisfied by definition, (16) follows immediately. □ 

So, we have shown the approximation property So, we obtain the 
following overall convergence result. 

Theorem 5 Assume that 


— the regularity assumptions (R) and (Rl) are satisfied on the domain 17, 

— the problem is discretized using the Taylor-Hood element and 

— the normal equation smoother introduced above is used as smoother. 

Then the two-grid method converges if sufficiently many smoothing steps 
are applied, i.e., we have 

- 4lllo.fc < 9 (j^)II| 4°^ - 4lllo.fc, 

with q{v) := CsCav~^^'^, where the constants Ca and Cs are independent of 
the grid level k and the choice of the parameter a. 
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Note that we have shown that the method converges in the norm ||| • |||o,fe if 
sufficiently many pre-smoothing steps are applied. The application of post¬ 
smoothing steps does not derogate the convergence because the proposed 
smoother is power-bounded. Moreover, if we assume that only post-smoothing 
steps are applied, the combination of smoothing property and approximation 
property (which now have to be combined in the inverse order) leads to con¬ 
vergence in the residual norm ||| • ||| 2 ,fe. Again, due to power-boundedness of 
the smoother, the method stays convergent if, besides sufficiently many post¬ 
smoothing steps, also pre-smoothing steps are applied. 

For all the mentioned cases, the convergence of the W-cycle multigrid 
method follows under weak assumptions, cf. [ 3 . 


5 Numerical Results 


In this section, we illustrate the convergence theory presented within this paper 
with numerical experiments. 

The domain 17 was chosen to be the unit square 17 := (0,1)^. As mentioned 
in Section]^ the weak inf-sup-condition (S) can be shown for the Taylor-Hood 
element only if at least one vertex of each element is located in the interior 
of the domain 17. As this is not satisfied for the standard decomposition of 
the unit square into two triangular elements, we have chosen the coarsest grid 
(grid level fc = 0) to be a decomposition of the domain 17 into 8 triangles, cf. 
Fig-El The grid levels k = 1,2,... were constructed by uniform refinement, 
i.e., every triangle was decomposed into four subtriangles. 

The desired velocity field (desired state) ud was chosen to be 






4 

5 


The desired velocity field is visualized in the left-hand-side picture in both. 
Fig. and 

For solving the discretized KKT-system, we have used the proposed W- 
cycle multigrid method. We have applied v pre- and v post-smoothing steps 
using the normal equation smoother. The matrix Ck was chosen as follows 


^ \ 
Sk 

a~^Ak 

V a~^SkJ 


(17) 


where A/. := diag {Mjj^k + a^^^ATc/.fc) and '■= adiag (DkAf^^Dl). Here, 
Mjj^k and Kjj^k are the mass matrix and the stiffness matrix, representing the 
L^-inner product and the i7^-inner product in Uk, respectively. The matrix 
Dk represents the bilinear form d{uk,Pk) = (V • Uk,Pk)L^(a) on Uk x Pk- Note 
that the matrix Ck, introduced above, is spectrally equivalent to the matrix Ck, 






18 


Stefan Takacs 


introduced in Section Therefore, the choice proposed above is also covered 
by the convergence theory. The damping parameter was chosen to be r = 0.35 
for all grid levels k and all choices of a. 



Fig. 1 Discretization on grid levels fc = 1 and /c = 2, where the squares denote the degrees 
of freedom of (the components of) u and A are the the dots denote the degrees of freedom 
of p and fi 
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Fig. 2 Desired velocity field ii£), optimal velocity field u and optimal control / for a = 1 
on grid level fc = 3 
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Fig. 3 Desired velocity field optimal velocity field u and optimal control / for a = 
10“^^ on grid level k = 3 


i 


The solution of the optimal control problem can be seen in Fig. and 
Note that the desired velocity field is a general L^-function (due to the 
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jump) and therefore it cannot be reached by the (optimal) velocity field which 
is an -function and therefore continuous. For the case ct = 1, we observe 
the optimal velocity held and the control to be rather smooth. (The control 
does not take large values). For small values of a, like a = 10“^^, the desired 
velocity held is approximated quite well, cf. Fig. This is achieved by rather 
large values of the control (high forces) which are concentrated on the region 
where the desired velocity held has its jump. Forces have to be applied with the 
same orientation as the desired state, as well as with the opposite orientation. 
(In the picture mainly the forces with opposite orientation can be seen.) As 
mentioned in the introduction, we are interested in a fast linear solver which 
also works well for such small choices of a. 

The number of iterations and the convergence rate were measured as fol¬ 
lows: we start with = 0 and measure the reduction of the error in each 
step using the residual norm ||| • ||| 2 ,fc. The iteration was stopped when the 
initial error was reduced by a factor of e = 10“®. The convergence rates q is 
the mean convergence rate in this iteration, i.e.. 


( - xlW2,k 

9 = ^ 

VIlFfc -2^fclll2.fe 


1 /n 


where n is the number of iterations needed to reach the stopping criterion. 

(i) 

Here, x*j, is the exact solution and x), is the z-th iterate. 


ly = 

■■1 + 1 

u = 

: 2 -1- 2 

ly = 

: 4-1-4 

ly = 

00 

+ 

00 

n 

Q 

n 

1 

n 

<? 

n 

<? 

61 

0.796 

32 

0.647 

21 

0.507 

15 

0.390 


Table 1 Number of iterations n and convergence rate q depending on = Upre + impost, the 
number of pre- and post-smoothing steps, on grid level k = 4 for a = 1 




a = 

n 

1 

Q 

a = 

n 

10-3 

q 

a. = 

n 

! 

to 1 

1 i 

s 

a = 

n 

10-9 

q 

OL = 

n 

10-12 

q 

k = 

3 

32 

0.648 

33 

0.651 

35 

0.673 

48 

0.749 

51 

0.760 

k = 

4 

32 

0.647 

32 

0.646 

33 

0.657 

46 

0.738 

73 

0.827 

k = 

5 

32 

0.645 

32 

0.644 

32 

0.644 

39 

0.697 

60 

0.793 

k = 

6 

31 

0.636 

31 

0.636 

31 

0.635 

32 

0.647 

46 

0.739 

k = 

7 

29 

0.620 

29 

0.620 

29 

0.618 

29 

0.621 

42 

0.716 


Table 2 Number of iterations n and convergence rate q for z/ = 2 -h 2 pre- and post¬ 
smoothing steps 
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In Table we compare for a fixed grid level (level /c = 4) and a fixed choice 
a = 1 the convergence rates for several choices of v, the number of pre- and 
post-smoothing steps. We see that the convergence rates behave approximately 
like . This is consistent with the theory which guarantees the convergence 
rate being bounded by C as this only describes the asymptotic behavior. 

In Tablej^we compare various grid levels k and choices of the parameter a. 
Here, we have used a fixed choice of i/ = 2 + 2 pre- and post-smoothing steps. 
First we observe that the number of iterations seems to be well-bounded for all 
grid levels k which yields an optimal convergence behavior. Moreover, we see 
that the number of iterations is also well-bounded for a wide range of choices 
of the parameter a, i.e., we observe also robust convergence as predicted by 
the convergence theory. 

It has to be mentioned that for the model problem, also the (more efficient) 
V-cycle multigrid method converges with rates comparable to the convergence 
rates of the W-cycle multigrid method. However, the V-cycle is not covered 
by the convergence theory. 


6 Conclusions and Further Work 

In the present paper we have shown that the construction of an all-at-once 
multigrid method for a Stokes control problem is possible. Here, a precondi¬ 
tioned normal equation smoother was chosen. The overall numerical complex¬ 
ity of this method seems to be comparable to block-preconditioned MINRES 
iterations, cf., e.g.. Table 4.2 in [18], which shows the number of MINRES 
iterations needed. (Note that in each MINRES step one multigrid cycle is ap¬ 
plied to each component of the overall block-matrix, i.e., to the velocity it, the 
pressure p, the adjoined velocity A and the adjoined pressure p,.) 

One advantage of the all-at-once multigrid method, introduced in the 
present paper, is the fact that an outer iteration is not necessary, the multigrid 
iteration is an linear iteration scheme which can be directly applied to solve the 
problem. As we could show the approximation property for a particular choice 
of norms, the construction of other smoothers is of particular interest. The 
convergence rates we have observed in this paper for a multigrid method with 
normal equation smoothing are comparable with the convergence rates ob¬ 
served in m for a multigrid method with normal equation smoothing applied 
to an optimal control problem with elliptic state equation. For that problem 
we have seen that other smoothers are available which lead to much faster 
convergence rates, cf. m and others. Similar improvements were possible for 
the generalized Stokes problem, cf. [14]. Therefore, it seems to be reasonable 
to construct faster smoothers also for the the Stokes control problem. 
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